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Abstract 

The coexistence between the hep and bcc phases of sohd ^He at fixed pressure, and at three different 
temperatures, is studied by means of the path integral Monte Carlo method. Microscopic results for the 
main energetic and structure properties of the two phases are reported. The differences between hep and 
bcc ^He are shown to be small with the obvious exception of the static structure factor. When crossing 
the phase transition line, the most significant changes are observed in the kinetic energy per particle and 
in the Lindemann ratio, both pointing to a less correlated quantum solid in the bcc side. No off-diagonal 
long-range order is observed in the hep and bcc perfect lattices. 
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I. INTRODUCTION 



Quantum solids have been studied for many years in the context of low-temperature condensed- 
matter physics. Contrarily to classical solids, the high zero-point motion of the atoms makes their 
oscillation around the lattice sites significant even at very low temperature. This makes quantum 
crystals have large Lindemann ratios and unambiguous anharmonic effects. Solid helium in its two 
stable isotopes ^He and ^He, and their mixtures, represent the best example of quantum solid [it]. 
Helium has a mass so small and its atomic interaction is so shallow that remains liquid even in 
the limit of zero temperature. It is necessary to pressurize ^He and ^He up to ~ 25 and 30 bar, 
respectively, to freeze them at T = 0. Because of the different mass and their different quantum 
statistics (^He is a boson and ^He is a fermion), the solidification is in the two cases different: ^He 
crystallizes into an hep solid whereas ^He does it into a bcc one. The phase diagram of both systems 
offers also the counterpart: ^He evolves to an hep phase at higher pressures and ^He presents a 
very narrow bcc region at temperatures T ~ 1.5 — 1.7 K. 

The theoretical and experimental interest in solid ^He has been reinforced since the experimental 
observation in 2004 of non classical rotation of inertia (NCRI) in torsional oscillator experiments 
by Kim and Chan [3]. This experiment, subsequent ones by the same team, and independent 
measures in a number of different laboratories have confirmed these findings but the fraction of mass 
decoupled (superfluid fraction) depends strongly on the quality of the crystal §]. The experimental 
data and their interpretation as unambiguous signatures of the existence of the supersolid state 
are still matter of debate. The observation by Beamish and collaborators [1] of an increase of the 
shear modulus of hep ^He at nearly the same onset temperature for NCRI (Tq ~ 100 mK) has lead 
to think that, at least in a partial way, the NCRI signal can be attributed to the change in the 
elastic properties. However, there is a general consensus that not all the torsional oscillator results 
can be interpreted solely on plasticity terms f5:, (jj. 

The main focus of this recent research activity has been hep ^He since this is the stable phase 
at very low temperatures. However, Eyal et al. [3] have reported recently torsional oscillator 
experiments on bcc ^He at temperatures between 1.3 and 1.9 K that show similar phenomena to 
those measured in the hep phase below 100 mK. This is a surprising result since, even if the disorder 
of the sample is eventually large, the temperatures of this experiment are one order of magnitude 
larger than in hep. This is not the only relevant result related to the bcc phase that emerged in 
the last years. Some time ago, Markovich et al. reported neutron-scattering measurements of 
phonons in bcc "^He and found an unexpected "opticlike" mode along the [110] direction. This new 
mode has been theoretically interpreted by Gov [9] in terms of correlated dipolar interactions. 

The boundaries of the bcc phase of ^He were first determined by Vignos and Fairbank [1(?] in 
1961 after analyzing their measurements of the velocity of longitudinal sound. The existence of this 
phase was confirmed after x-ray scattering by Schuch and Mills 11]. Later on, Grilly and Mills [l2|l 
measured the volume change along the hcp-bcc coexistence line and Edwards and Pandorf [13] 
performed accurate measurements of the heat capacity. Updated accounts of the thermodynamic 
properties of bcc ^He can be found in Refs. [id . [3]. In Fig. [H the P-T phase diagram in the 
region of interest is plotted. As one can see, the stability of the bcc solid is only ~ 0.04 K wide at 
fixed pressure and it represents an intermediate phase between the hep solid and the superfluid in 
the range of temperatures T ~ 1.5-1.75 K. 

From the theoretical point of view, the high anharmonicity of solid ^He makes its analysis 
much more demanding since the usual harmonic theory used in classical solids is not applicable. 
The self-consistent phonon (SCP) theory incorporates approximately the short-range correlations 
induced by the interatomic potential, leading to reasonable descriptions of many of its properties [l|. 
However, a more accurate study of quantum solids requires of a microscopic approach. To this end, 
quantum Monte Carlo (QMC) methods have been widely applied in the past to the study of solid 
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FIG. 1: (Color online) Pressure-temperature phase diagram in the stable region of the bcc phase of solid 
^He. The squares correspond to the thermodynamic points studied in the present work. 



^He. However, the main part of the QMC work is devoted to the hep phase since it is the stable 
crystal in the limit of very low temperature. Results on the bcc phase are scarce and concentrated 



to the limit of zero temperature |l6l4l8l] . To our knowled ge, the only existing work in the bcc 
stable P-T region shown in Fig. [1] was done by Ceperley [l9| with the main goal of estimating 
the momentum distribution and kinetic energies after some neutron scattering experiments (i^. 
Our aim in the present work is to carry out a QMC study of the bcc phase at three temperatures 
and densities corresponding to the experimental points shown in Fig. [T] and compare them with 
the corresponding hep ones along the coexistence line. As we will show, the differences between 
both phases are found small except for the kinetic energies and Lindemann ratios where the less 
correlated character of the bcc lattice is revealed. 

The rest of the paper is organized as follows. In Sec. II, we introduce the path integral Monte 
Carlo (PIMC) method used in our simulations. Sec. Ill comprises the main results obtained, and 
finally a brief account of the main conclusions is contained in Sec. IV. 



II. THEORETICAL APPROACH 



The structure and energetic properties of the hcp-bcc coexistence line have been studied using 
the path integral Monte Carlo (PIMC) method. PIMC provides a fundamental approach in the 
study of strongly interacting quantum systems at finite temperature by means of a stochastic 
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estimation of the thermal density matrix [2]}]. As it is well known, the partition function 

Z = TV(e-^^) = y" dR(R|e-^^|R) (1) 

allows for a full microscopic description of the properties of a given system with Hamiltonian 
i7 = + at a temperature T = {ksji)'^ (we use the position basis \R) = |ri, . . . ttv), with 
N the number of particles). However, in quantum systems the noncommutativity of the kinetic 
and potential energy operators (respectively, K and V) makes impractical a direct calculation of 
Z using its definition ([T]). 

The basic idea of PIMC is to use the convolution property of the thermal density matrix 
p{R, R'; P) = {R\e~^^\R'), in order to rewrite the partition function as 

Z = / n dR,p{Ri,R,+i-£) , (2) 

i=0 

with e = P/M and the boundary condition Rm = Ro- For sufficiently large M, we recover 
the high-temperature limit of the thermal density matrix, where it is legitimate to separate the 
kinetic contribution from the potential one (Primitive Approximation). If one ignores the quantum 
statistics of particles, the distribution law appearing in Eq. [2] is positive definite and can be 
interpreted as a probability distribution function which can be sampled by standard metropolis 
Monte Carlo methods. 

In practice, the PIMC method consists in mapping the finite-temperature quantum system 
onto a classical system made up of closed ring polymers. This technique may be referred as 
an "exact" method, in the sense that using an accurate approximation for the high-temperature 
density matrices, the results are not affected by this approximation within the statistical error. 
However, the number M of convolution terms (beads) necessary to reach the convergence of Eq. 
[2] to the exact value of Z is inversely proportional to the temperature of the system. This means 
that, when approaching the interesting quantum regime at very low temperature, M increases fast 
making simulations hard, if not impossible, due to the very low efficiency in the sampling of the long 
chains involved. To overcome this problem, it is important to work out high-order approximation 
schemes for the density matrix, able to work with larger values of e. The approximation we use in 
this work is called Chin Approximation (CA) |22|. CA is based on a fourth-order expansion of the 
operator e~^^ which makes use of the double commutator [[V, K],V], this term being related to the 



gradient of the interatomic potential. With respect to the Takahashi-Imada approximation [23|, 
which is accurate to fourth order only for the trace, the new feature appearing in the CA is 
the presence of coefficients weighting the different terms in the expansion of the action. These 
coefficients are continuously tunable, making possible to force the error terms of fourth order to 
roughly cancel each other and get an effective sixth-order approximation. 

An additional problem one has to deal with when simulating quantum many-body systems with 
PIMC arises from the indistinguishable nature of particles. In the case of bosons like ^He, the 
indistinguishability of particles does not change the positivity of the probability distribution in Eq. 
[2] and the symmetry of p{R, R'; f3) can be recovered via the direct sampling of permutations between 
the ring polymers representing the quantum particles. To this purpose, we use the Worm Algorithm 
(WA) ^] which samples very efficiently the permutation space. The basic idea of this technique is 
to work in an extended configuration space, given by the union of the ensemble Z, formed by the 
usual closed-ring configurations, and the ensemble G, which is made up of configurations where all 
the polymers but one are closed. Thanks to the presence of an open polymer, one can search the 
atoms involved in a permutation cycle by means of single particle updates, which do not suffer of a 
low acceptance rate and guarantee an efficient and ergodic sampling of the bosonic permutations. 
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TABLE I: PIMC results of the total {E/N) and partial energies [V/N, K/N) for the three thermodynamic 
points in the hcp-bcc coexistence line here analyzed. Experimental values of the bcc energies, reported in 
the first row, are taken from Ref. Figures in parenthesis stand for the statistical errors. 
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We have to notice that the probability distribution used to sample the configurations in G is not 
equal to the one appearing in Eq. [2] and, therefore, these configurations cannot be used to calculate 
diagonal properties, such as the energy or the superfluid density. However, the G-configurations 
can be used to compute ofF-diagonal observables such as the one-body density matrix /9i(ri, r'j^). 
Furthermore, the WA samples both diagonal and off-diagonal configurations and then it is able 
to give an estimation of the normalization factor of /9i(ri,r^). In this way, one can compute the 
properly normalized one-body density matrix and so avoid the systematic uncertainties introduced 
by an a posteriori normalization factor. 



III. RESULTS 

We have carried out our PIMC study of the hcp-bcc coexistence line particularizing the simula- 
tions in the three experimental thermodynamic points shown in Fig. [H Explicitly, in the P-T phase 
diagram the coordinates of these points are A= (1.5, 26.419), B= (1.6, 27.572), and C= (1.7, 29.029) 
with units of bar and kelvin for the pressure and temperature, respectively. The corresponding 
densities are taken from the experimental data contained in Refs. 

B, Q: {p'r = 0.028834, 

pY = 0.028571), {p'g'P = 0.028954, p)^"" = 0.028679), and {p^^^ = 0.029080, /oj^f'^ = 0.028805), all 
in units of . 

In the Hamiltonian of the system the potential part is built as a sum of pair-wise interatomic 



potentials, V = J2i<j ^(^ij)) with V{r) of Aziz type [25|. Simulations are performed using periodic 



boundary conditions with a number of atoms per simulation cell of N = 180 (hep) and N = 128 
(bcc); tail corrections to the energy due to the use of a finite number of particles are estimated by 
running simulations with different values and extrapolating its linear behavior in 1/A^ to the 
thermodynamic limit 1/A^ ^ [2^. The convergence to the limit e — t- is achieved with relatively 
large values of e (e = 0.010 K~^) due to the high accuracy of the CA [22]. 

PIMC results for the total energy per particle in the three mentioned thermodynamic points 
are reported in Table 1. The difference between the energies of the two lattices in the coexistence 
line is small but, in the three cases, the bcc energies are larger than the hep ones (in absolute 
values). This difference is mainly due to the slightly smaller density in the bcc side with respect 
to the hep one, with a smaller effect of the lattice type (at the same density). On the other hand, 
PIMC energies lie below the available experimental energies of the bcc solid in an amount ~ 0.20 
K, a feature already observed in liquid ^He and attributable to the particular Aziz potential used 



in the present simulations 27| 



In Table 1, we also have included the partial energies per particle. The results for the kinetic 
energy per particle K/N are the most relevant because it is possible to measure them using deep 
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FIG. 2: (Color online) Density profile /i(r) of ^He atoms around the lattice sites for the two lattices at 
T = 1.5 K. 



inelastic neutron scattering. As it is well known, if the momentum of the incoming neutron is high 
enough the impulse approximation holds and the atomic momentum distribution n(/c) is attainable, 
and from it the kinetic energy. Blasden et al. reported experimental results ol K/N for both 
solid lattices: at T = 1.07 K and hep K/N = 23.6 K, and at T = 1.72 K and bcc K/N = 23.7 K 
(the densities of both phases were the same). If just the classical estimation 3/2kT is considered for 
an isothermical comparison, one can see that the kinetic energy of the bcc crystal is significantly 
smaller. PIMC results by Ceperley [19| also agree with this trend (in this case the comparison is 
made between fee and bcc lattices). Our present results allow for a best comparison since they 
are obtained following the coexistence line. As one can see in Table 1, there is a nearly constant 
decrease of ~ 0.5 K when going from the hep point to the bcc one and our results are in overall 
agreement with the experimental determinations of Blasdell et al. |20l |. The thermal effects in the 
kinetic energy for any of the two lattices are dominated by the classical term 2>/2kT that gives a 
difference of 0.15 K from one point to the next one. 

We have also calculated some structural properties of the bcc phase to better characterize it and 
to identify specific signatures of this lattice symmetry. One of the most relevant in the microscopic 
study of any solid phase is the density profile /u(r), defined as the probability of finding a particle 
in the interval (r, r + dr) around any of the lattice points of the crystal. Results of //(r) at the 
thermodynamic point A and for the two coexistence lattices are shown in Fig. [2j As expected, the 
density profile of the bcc crystal is a bit wider than the hep one, with a decrease of localization as it 
corresponds to its slightly more open structure. From the density profile it is possible to calculate 
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FIG. 3: (Color online) Two-body radial distribution function for the hep and bee phases of solid "^He at 
T = 1.5 K. 



the mean squared displacement around a site, 



47r 



dr r^fi{r) 



and from it to estimate the Lindemann ratio, 
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(3) 



(4) 



a being the nearest-neighbor distance in the perfect crystalline lattice. We have obtained 7hcp = 
0.26 and 7bcc = 0.28, in agreement with experimental data [ij. The Lindemann ratio is a good 
measure of the zero-point motion and thus the larger value of 7 for the bcc lattice points to an 
enhancement of its quantum nature. 

In Fig. [3l we show results of the two-body radial distribution function for both lattices in the 
coexistence at T = 1.5 K. One can see the usual oscillations due to the periodic spatial order and 
only a tiny difference between both lattices. The height of the main peak of the bcc g{r) is slightly 
smaller than the hep one but this is more an effect of its smaller density than a consequence of 
the lattice symmetry. The difference in the spatial order induced by the lattice is more clear in 
reciprocal space by defining the static structure factor 




q (A ) 



FIG. 4: (Color online) Static structure factor for the hep and bee phases of solid "^He at T = 1.5 K. 



Results for S{q) at both sides of the coexistence line at T = 1.5 K are shown in Fig. HI Several 
Bragg peaks are clearly identified; they are located at the expected points for each type of lattice. 
As it is well known from solid state physics, the determination of the static structure factor is the 
best way for identifying the particular lattice of the crystal and, in fact, the x-ray measures by 
Schuch and Mills 1]| were the definite prove of the existence of a stability region for bcc ^He. 



Much of our present interest in the bcc phase of solid ^He lies on the unexpected results found 
by Eyal et al. [l\ which point to disorder-induced mobility, that manifests in the TO experiments 
as the typical mass decoupling observed in hep at very low temperatures (T < 200 mK). This could 
open the possibility of a supersolid scenario at higher temperatures. Any signature of supersolidity 
would appear in off-diagonal properties of the system, like the one-body density matrix pi{r) 
defined as 



Pi(ri,r'i) 



1 

Z 



dT2, ...,rNp{r 



1, 



,rAr;ri, 



>rAr;/3) . 



(6) 



With the help of the WA one can calculate this function with high efficiency and avoid normalization 
problems that were present in other permutation sampling methods. In Fig. [5l we show results for 
pi(r) at T = 1.5 K calculated for the two lattices. The shape of the two results is very similar and 
coincides with previous estimations: a fast decay up to ~ 4 A, a small increase near the position 
of the first neighbor, and finally a kind of exponential decay to zero. The most relevant feature 
of these results is the absence of off-diagonal long-range order, which would manifest as a nonzero 
asymptotic constant (condensate fraction) when r oo. Therefore, neither the perfect hep crystal 
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nor the perfect bcc one show signals of supersoHdity according to our PIMC results. 



IV. CONCLUSIONS 

We have performed a microscopic study of the coexistence hcp-bcc line in solid ^He using the 
path integral Monte Carlo method. This method allows for a nearly ab-initio study where the only 
inputs are the mass and the interatomic potentials. Our goal has been to get an accurate view of 
the possible differences between both lattices along the measured coexistence line. Apart from the 
intrinsic interest on the study of the bcc phase, which has been scarcely studied in the past, we 
were stimulated by the recent activity of Polturak and collaborators [3] who have shown intriguing 
results using the TO technique. 

Our results show a small influence of the lattice type when crossing the coexistence line. The 
internal energies are similar and the differences between both crystals are more attributable to the 
slight difference in density than to the particular crystal symmetry. Also the radial distribution 
functions show only very tiny differences. What is more significant is the drop in the kinetic 
energy from the hep phase to the bcc one, that we have estimated to be 0.5 K, and that is in 
overall agreement with experimental findings [20]. The less compact structure of the bcc solid is 
also observed in the results obtained for the density profiles and the Lindemann ratios. Finally, we 
have calculated the one-body density matrix in both sides of the coexistence line without observing 
any relevant difference between them and, more importantly, without obtaining off-diagonal long- 
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range order which allowed us to conclude that a perfect bcc crystal, at the temperatures where this 
phase is stable, is not a supersolid. Work is in progress in our group to determine if the introduction 
of point defects (vacancies) in the bcc crystal could induce supersolidity at these apparently too 
high temperatures. 
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